#!/bin/bash

export LC_ALL=C

export BASE_DATA_PATH=/home/mixagol/data/
export BASE_DATE=20120802
export CUR_DIR=${BASE_DATA_PATH}/5_pathways/${BASE_DATE}

mkdir -p $CUR_DIR

##
## Получаем описание всех путей
##
5_pathways/get_full_pathways.py > ${CUR_DIR}/_pathways_full.txt
sort -t$'\t' -k1,1 -o ${CUR_DIR}/_pathways_full.txt ${CUR_DIR}/_pathways_full.txt

##
## Получаем содержание всех путей
##
mkdir ${CUR_DIR}/tmp_input/
mkdir ${CUR_DIR}/tmp_output/
cat ${CUR_DIR}/_pathways_full.txt | cut -f1 | split -l 1000 - ${CUR_DIR}/tmp_input/

export PROCESSES_NUM=10
ls ${CUR_DIR}/tmp_input/ \
    | xargs -n1 -P $PROCESSES_NUM -I {} bash -c \
    "cat ${CUR_DIR}/tmp_input/{} | 5_pathways/get_full_kegg_pathway_map.py > ${CUR_DIR}/tmp_output/{}"

cat ${CUR_DIR}/tmp_output/* | sort -t$'\t' -k1,1 > ${CUR_DIR}/pathway_kegg_map_full.txt
cat ${CUR_DIR}/pathway_kegg_map_full.txt | awk -F'\t' '{print $2"\t"$1}' | sort -t$'\t' -k1,1 > ${CUR_DIR}/kegg_pathway_map_full.txt

rm -rf ${CUR_DIR}/tmp_input/
rm -rf ${CUR_DIR}/tmp_output/

##
## Добавляем число генов в каждом пути
##
join -t$'\t' \
    <(cat ${CUR_DIR}/pathway_kegg_map_full.txt | cut -f1 | uniq -c | awk '{print $2"\t"$1}') \
    ${CUR_DIR}/_pathways_full.txt \
    > ${CUR_DIR}/pathways_full.txt

rm -rf ${CUR_DIR}/_pathways_full.txt 

##
## Получаем соответствие наших ncbi_id -> kegg_id
##
mkdir ${CUR_DIR}/tmp_input/
mkdir ${CUR_DIR}/tmp_output/
cat ${BASE_DATA_PATH}/1_databases/${BASE_DATE}/genes_info.txt | cut -f2 | split -l 10000 - ${CUR_DIR}/tmp_input/

export PROCESSES_NUM=10
ls ${CUR_DIR}/tmp_input/ \
    | xargs -n1 -P $PROCESSES_NUM -I {} bash -c \
    "cat ${CUR_DIR}/tmp_input/{} | 5_pathways/create_ncbi_kegg_map.py > ${CUR_DIR}/tmp_output/{}"

cat ${CUR_DIR}/tmp_output/* | sort -t$'\t' -k1,1 > ${CUR_DIR}/ncbi_kegg_map.txt

rm -rf ${CUR_DIR}/tmp_input/
rm -rf ${CUR_DIR}/tmp_output/

##
## Получаем kegg_id -> gene definition
##
mkdir ${CUR_DIR}/tmp_input/
mkdir ${CUR_DIR}/tmp_output/
cat ${CUR_DIR}/ncbi_kegg_map.txt | cut -f2 | sort | split -l 10000 - ${CUR_DIR}/tmp_input/

export PROCESSES_NUM=10
ls ${CUR_DIR}/tmp_input/ \
    | xargs -n1 -P $PROCESSES_NUM -I {} bash -c \
    "cat ${CUR_DIR}/tmp_input/{} | 5_pathways/get_genes_info.py > ${CUR_DIR}/tmp_output/{}"

cat ${CUR_DIR}/tmp_output/* | sort -t$'\t' -k1,1 > ${CUR_DIR}/kegg_genes_functions.txt

rm -rf ${CUR_DIR}/tmp_input/
rm -rf ${CUR_DIR}/tmp_output/

join -t$'\t' -1 2 -2 1 \
    <(cat ${CUR_DIR}/gene_kegg_map.txt | cut -f1,2 | sort -t$'\t' -k2,2) \
    <(cat ${CUR_DIR}/kegg_genes_functions.txt) \
    | cut -f2,3 \
    | sort -t$'\t' -k1,1 \
    > ${CUR_DIR}/gene_functions.txt



##
## Строим соотвтствие gene_id -> ncbi_id
##
join -t$'\t' -1 2 -2 1 \
        <(cat ${BASE_DATA_PATH}/1_databases/${BASE_DATE}/genes_info.txt | cut -f1,2 | sort -t$'\t' -k2,2) \
        <(cat ${CUR_DIR}/ncbi_kegg_map.txt | cut -f1,2 ) \
    | cut -f2,3 \
    | sort -t$'\t' -k1,1 \
    > ${CUR_DIR}/gene_kegg_map.txt

##
## Строим соотвтствие gene_id -> pathway
##
join -t$'\t' -1 2 -2 1 \
    <(cat ${CUR_DIR}/gene_kegg_map.txt | sort -t$'\t' -k2,2 ) \
    ${CUR_DIR}/kegg_pathway_map_full.txt \
    | cut -f2,3 \
    | sort -t$'\t' -k1,1 \
    > ${CUR_DIR}/gene_pathway_map.txt


##
## Строим pathway -> gene_vector
## Строим pathway -> avg_gene_vector
##
for E_VALUE in 1e+01 1e+00 1e-01 5e-02 1e-02 1e-03 1e-04 1e-05 1e-06 1e-07 1e-08 1e-09 1e-10 1e-11 1e-12 1e-13 1e-14 1e-15 1e-20 1e-30 1e-40 1e-50; do
    join -t$'\t' -1 2 \
        <(join -t$'\t' -2 2 \
            ${CUR_DIR}/gene_pathway_map.txt   \
            ${BASE_DATA_PATH}/3_normed_matrix/map_int_gene.txt  \
            | cut -f2,3 | sort -t$'\t' -k2,2 )\
        <(cat ${BASE_DATA_PATH}/3_normed_matrix/matrix_${E_VALUE}_.txt \
            | sort -t$'\t' -k1,1 -S200M -T ${CUR_DIR} )\
        | cut -f2,3 | sort -t$'\t' -k1,1 -S200M -T ${CUR_DIR} \
        | gzip \
        > ${CUR_DIR}/pathway_vectors_${E_VALUE}_.txt.gz
    
    zcat ${CUR_DIR}/pathway_vectors_${E_VALUE}_.txt.gz \
        | 5_pathways/create_avg_vectors.py \
        > ${CUR_DIR}/pathway_vectors_avg_${E_VALUE}_.txt
done

for Z_VALUE in 10 20; do
    join -t$'\t' -1 2 \
        <(join -t$'\t' -2 2 \
            ${CUR_DIR}/gene_pathway_map.txt   \
            ${BASE_DATA_PATH}/3_normed_matrix/map_int_gene.txt  \
            | cut -f2,3 | sort -t$'\t' -k2,2 )\
        <(cat ${BASE_DATA_PATH}/3_normed_matrix/matrix_${Z_VALUE}_.txt \
            | sort -t$'\t' -k1,1 -S200M -T ${CUR_DIR} )\
        | cut -f2,3 | sort -t$'\t' -k1,1 -S200M -T .\
        | gzip \
        > ${CUR_DIR}/pathway_vectors_wunw_${Z_VALUE}_.txt.gz
    
    zcat ${CUR_DIR}/pathway_vectors_wunw_${Z_VALUE}_.txt.gz \
        | 5_pathways/create_avg_vectors.py \
        > ${CUR_DIR}/pathway_vectors_avg_wunw_${Z_VALUE}_.txt
done

